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ABSTRACT 


The  objective  of  this  work  is  to  help  improve  seismic  monitoring  technology  through  the  development  and 
application  of  advanced  multivariate  inversion  techniques  to  generate  realistic,  comprehensive,  and  high-resolution 
3D  models  of  the  seismic  structure  of  the  crust  and  upper  mantle  that  satisfy  independent  geophysical  datasets.  Our 
focus  is  on  the  region  surrounding  Iran  from  the  cast  coast  of  the  Mediterranean  in  the  west,  to  Pakistan  in  the  east, 
an  area  of  prime  importance  to  NEM,  and  a  region  with  adequate  calibration  events  to  validate  our  model  and  to 
quantify  its  accuracy.  Specifically,  we  arc  working  to  integrate  surface-wave  dispersion,  receiver  function,  and 
satellite  and  ground-based  gravity  observations  to  help  constrain  the  shallow  seismic  structure  in  the  Arabian- 
Eurasian  collision  zone.  Building  on  our  earlier  work  combining  receiver  functions  and  surface  wave  dispersion,  and 
surface-wave  dispersion  and  gravity,  we  plan  to  continue  to  integrate  geophysical  data  sets  to  create  more  broadly 
compatible  earth  models.  We  also  explore  geologically  based  smoothness  constraints  to  help  resolve  sharp  features 
in  the  underlying  shallow  3D  structure. 


20101025206 


2010  Monitoring  Research  Review:  Ground-Based  Nuclear  Explosion  Monitoring  Technologies 


OBJECTIVES 

The  National  Nuclear  Security  Administration  (NNSA)  and  Air  Force  Research  Laboratory  (AFRL)  have  decided  to 
investigate  3D  modeling  as  part  of  the  effort  to  improve  knowledge  of  Earth's  congressional  and  shear  velocity 
structure.  Such  knowledge  will  help  reduce  uncertainty  in  our  ability  to  accurately  detect,  locate,  and  identify  small 
(mb  <  3.5)  seismic  events,  which  in  turn  will  lead  to  improved  capabilities  for  nuclear  explosion  monitoring  (NEM). 
For  seismieally  active  areas,  with  good  ground- truth  event  coverage,  earth  models  with  limited  accuracy  ean  be 
corrected  by  interpolating  results  from  nearby  ‘ground  truth'  events  (using  the  kriging  methodology)  making  it 
possible  to  detect,  locate,  and  identify  events  even  with  limited  resolution  of  Earth's  structure.  However,  sueh 
approaches  arc  less  effective  for  smaller  events,  and  event  location  and  characterization  remains  a  challenge  for 
aseismie  regions.  To  improve  monitoring  capability  in  sueh  instances,  we  must  develop  better  seismic  earth  models. 

The  objective  of  this  work  is  to  help  improve  seismic  monitoring  technology  through  the  development  and 
application  of  advanced  multivariate  inversion  techniques  to  generate  realistic,  comprehensive,  and  high-resolution 
3D  models  of  the  seismic  structure  of  the  crust  and  upper  mantle  that  satisfy  independent  geophysical  datasets.  Our 
focus  is  on  the  region  including  and  surrounding  Iran  (Figure  1 )  from  the  east  coast  of  the  Mediterranean  in  the  w  est, 
to  Pakistan  in  the  east,  an  area  of  prime  importance  to  NEM,  and  a  region  with  adequate  calibration  events  to 
validate  our  models  and  to  quantify  their  accuracy. 

Background 

Estimating  subsurface  geologie  variations  is  a 
challenge.  Seismologists  have  worked  on  the  problem 
for  more  than  a  century  (e.g ,  Milne,  1899;  Macelwane 
and  Sohon,  1936;  Dahlen  and  Tromp,  1999).  As  data 
quantity  has  increased  and  data  quality  and 
computational  ability  have  improved,  we  have  made 
important  advances  in  our  understanding  of  the 
subsurface.  Our  best  knowledge  applies  to  the 
shallowest  regions  as  well  as  depths  with  the  sharpest 
global  interfaces  (sediment-basement  contacts,  the 
base  of  the  crust,  base  of  the  mantle,  and  transitions 
near  410  km  and  660  km  depths),  where  resolution  is 
improved  as  a  result  of  the  strong  interactions  of 
seismic  body  waves  with  sharp  geologie  transitions 
(e.g.,  Helmberger,  1968;  Langston,  1979;  Shearer, 

1991 ;  Lay  et  al.,  2004).  We  have  also  done  well 
modeling  regions  with  smooth  velocity  changes  sueh 
as  the  lower  mantle,  which  allows  us  to  exploit  the 
information  in  teleseismie  body-wave  travel  times  to 
locate  seismic  sources  reasonably  well  (eg.,  Kennett, 
1991).  Still,  many  details  within  and  just  beneath  the 
lithosphere  elude  us.  We  have  been  able  to  surmise 
that  geologic  variations  here  are  substantial,  and  we 
know  that  they  frustrate  attempts  to  use  robust 
observations  sueh  as  regional  seismic  travel  times  to 
locate  events  in  many  parts  of  the  Earth  (e  g.,  Bondar 
et  al.,  2004). 

Travel-time  based  tomography  (e.g.,  Nolet,  2008)  opened  the  doors  to  3D  imaging  but  the  models  remain  blurry, 
often  suffer  from  interpretational  ambiguity,  and  are  not  easily  used  to  predict  other,  independent  seismic 
observations.  From  our  own  analyses  (Maeeira  et  al.,  2005;  Maeeira  and  Ammon,  2009),  we  have  seen  how  high- 
resolution  surface-wave  tomography  fails  to  produce  the  extremes  in  seismic  shear-wave  speed  that  arc  evident  from 
independent  observations.  In  particular,  achieving  a  model  w  ith  low  enough  seismic  wave  speeds  within  the  Tarim 
Basin  to  match  seismograms  from  high-quality  observations  remains  an  issue.  Waveform  tomography  methods 
improve  the  situation  somewhat,  including  information  from  both  the  amplitude  and  phase  of  the  signal,  but 
restriction  of  these  methods  to  lower  frequency  bands  limits  the  resolution  of  the  methods  and  the  substantial 


Figure  1  Map  of  focus  region  show  with  topographic 
and  bathymetric  shading  and  moderate  to  large 
earthquake  locations  (magnitudes  >  3.5  from  1990  to 
Spring  2008).  The  region  contains  the  Arabian  plate 
and  the  middle  segment  of  the  Alpine  to  Himalyan 
collision  zone,  which  is  constructed  primarily  of 
Phanerozoie  terranes  amalgamated  onto  southern 
Eurasia  during  the  closing  of  the  Tcthys  Ocean. 
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computation  requirements  of  these  methods  limit  their  application.  More  recent  finite- frequency  methods  (e.g.,  Zhou 
et  al,,  2004;  Dahlen  and  Zhou,  2006,  Nolet,  2008)  and  adjoint  methods  (e.g.,  Tarantola,  1984;  Tromp  et  al.,  2008) 
offer  more  complete  approaches  to  computing  sensitivity  kernels.  But  even  these  approaches  faee  limits  imposed  by 
data  bandwidth.  In  any  event,  such  fully  3D  waveform  methods  eould  benefit  greatly  from  aeeurate,  if  approximate, 
starting  models  derived  from  more  piecewise  interpretation  of  seismic  observables  combined  with  other 
observations. 

RESEARCH  ACCOMPLISHED 

Much  of  our  effort  during  the  first  year  of  this  project  has  been  on  the  software  development,  although  we  have 
completed  preliminary  inversions  that  we  present  below  with  little  discussion.  At  this  point  we  eaution  readers  about 
using  or  interpreting  these  preliminary  earth  models  they  will  improve.  We  begin  with  a  simple  conceptual 
illustration  of  the  challenges  we  address  and  eonelude  with  a  discussion  of  accomplishments  and  plans  for  the 
seeond  year  of  work. 

Combining  Gravity,  Rayleigh- Wave  Dispersion,  and  P-vvave  Receiver  Function  Observations 

Inversion  of  surface  wave  dispersion  observations  is  a  standard  method  for  estimating  3D  shear  velocity  structure  of 
Earth’s  erust  and  upper  mantle.  Nevertheless,  it  is  well  known  that  traditional  state-of-the-art  inversion  techniques 
suffer  from  poor  resolution  and  nonumqueness,  especially  when  a  single  surface-wave  mode  is  used  (Huang  et  al., 
2003).  This  is  particularly  true  at  shallow  depths  where  the  shorter  periods,  whieh  are  primarily  sensitive  to  upper 
crustal  structures,  are  difficult  to  measure  and  especially  true  in  teetonieally  and  geologically  complex  areas.  On  the 
other  hand,  regional  gravity  inversions  have  the  greatest  resolving  power  at  shallow  depths  sinee  gravity  anomalies 
decrease  in  amplitude  and  wavenumber  bandwidth  with  increasing  depth.  Gravity  measurements  also  supply 
constraints  on  rock  density  variations.  Thus  by  combining  surfaee-wave  dispersion  and  gravity  observations  in  a 
single  inversion,  we  can  obtain  a  self-eonsistent  high-resolution  3D  shear-veloeity/density  model  with  increased 
resolution  of  shallow  geologic  structures.  Receiver  function  analysis  (Langston,  1979)  is  a  well-established  tool  for 
imaging  relatively  sharp  changes  in  the  subsurface  structure  that  nicely  complements  information  contained  in 
surfaee-wave  dispersion  (e.g.,  Julia  et  al.,  2000). 

To  improve  our  view  of  the  earth  structure  within  the 
focus  region,  we  are  working  to  eombine  receiver 
functions  with  surfaee-wave  dispersion  and  gravity 
observations.  The  combination  of  receiver  functions 
with  tomographieally  localized  surfaee-wave 
dispersion  is  well  established  (e.g.,  Julia  et  al.  2000) 
The  integration  of  gravity  observations  with 
surfaee-wave  dispersion  is  a  more  recent  development 
(e.g.,  Maeeira  and  Ammon,  2009),  but  these  data  are  a 
good  match.  To  a  large  part,  the  gravity  observations 

add  important  constraints  on  the  location  of  strong 
shallow  heterogeneity  to  the  information  on  smooth 
variations  found  in  period-dependent  dispersion  maps. 
The  choiee  of  a  density- velocity  relationship  ean  be 
subjective  (Julia  et  al.,  2004;  Broehcr,  2005;  Maeeira 
and  Ammon,  2009),  but  testing  the  sensitivity  of  the 
results  to  a  range  of  relations  is  not  a  challenge.  To 
isolate  gravity  signatures  associated  with  density 
variations  (as  opposed  to  dynanne  and  flexure 
variations)  we  plan  to  filter  the  gravity  signals  to 
inelude  short  wavenumbers  sensitive  to  the  erust. 

These  concepts  arc  summarized  conceptually  in 
Figure  2,  whieh  shows  the  regions  of  the  lithosphere 
most  sensitive  to  the  different  data  that  we  employ. 

To  eonstruet  an  approximate  3D  model  of  the  lithosphere,  we  use  a  hybrid  I D-3D  inversion.  In  many  tomography 
analyses,  dispersion  variations  are  converted  to  shear-veloeity  variations  by  inverting  dispersion  curves  extracted 
from  the  tomographic  model  for  a  localized  1 D  structure.  Smoothness  constraints  are  applied  to  reduce  eell-to-eell 


Figure  2  Conceptual  illustration  of  the  primary 
sensitivities  of  different  data  sets  to  the  different  depths 
of  the  lithosphere  in  a  simultaneous  inversion  of 
surfaee-wave  dispersion,  gravity,  and  receiver-function 
observations.  The  crust  shown  contains  a  large  basin 
(yellow).  The  receiver  functions  provide  information  on 
strong  velocity  contrasts  such  as  the  erust-mantle 
transition  and  near-surfaec  structures;  the  spatially 
filtered  grav  ity  provides  constraints  on  near-surfaee 
structures;  and  dispersion  provides  information  on  the 
absolute  seismic  velocities  throughout  the  structure. 
Note  that  this  cartoon  is  conceptual,  the  sensitivity  of 
the  data  is  more  subtle  than  shown. 
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shear-velocity  fluctuations  in  the  resulting  composite  3D  structure.  Since  gravity  observations  can  be  efficiently 
modeled  using  prisms  they  provide  direct  information  on  eell-to-eell  density  variations.  In  a  particular  cell,  surface 
wave  dispersion  can  be  inverted  simultaneously  with  receiver  functions  (surface  or  downward-continued  signals)  In 
regions  with  receiver  function  overlap,  multiple  signals  can  be  combined  to  produce  a  structure  that  matches 
compatible  features  in  several  receiver  functions.  Because  of  the  hybrid  nature  of  the  inversion,  the  resulting  model 
is  an  approximation  to  the  true  3D  structure,  but  we  can  use  to  provide  a  starting  point  for  future  full  3D 
waveform-based  inversions  of  the  subsurface. 

Preliminary  Application  to  the  Iranian  Plateau  and  Surrounding  Regions 


2  km  10  km  20  km 


Figure  3.  Results  of  a  preliminary  inversion  of  Rayleigh-wave  dispersion  and  Bouguer  gravity 
observations  to  estimate  shear-veloeity  beneath  the  foeus  region.  The  inversion  was 
performed  using  the  existing  SVD-based  approach  and  took  about  36  hours  to  compute  on  a 
high-end  workstation.  This  result  is  preliminary  and  should  not  be  used  for  specific 
interpretation.  The  extreme  slow  values  for  the  shallow  Caspian  Sea  are  obvious;  deeper 
values,  particularly  those  with  high  wavenumbers  are  likely  artifacts.  The  inversion  weighted 
observations  at  shorter  periods  much  more  than  longer-period  observations. 

A  large  component  of  the  proposed  work  involves  improving  current  software  to  allow  more  easy  extension  to  larger 
foeus  areas  and  more  flexible  incorporation  of  additional  data.  For  example,  the  earlier  work  in  Maeeira  and  Ammon 
(2009)  blended  gravity  and  surfaec-wave  observations  to  produce  a  shear-veloeity/density  model  for  central  Asia, 
but  did  not  include  receiver-function  information.  Our  inversion  is  a  straightforward  composite  of  coupled  and 
uncoupled  linearized  inversions  of  nonlinear  data  and  model  relationships.  The  preliminary  results  we  present  below 
do  not  have  the  gravity  wavenumber  filters  applied.  We  rely  on  the  dispersion  models  from  the  University  of 
Colorado  group  (Ritzwoller  and  Levshin.  1998;  Levshin  et  al.,  2001,  2002)  for  the  preliminary  results,  but  we  will 
include  other  measures  as  the  project  proceeds.  Figure  3  is  a  plot  of  first-order  shear-velocity  variations  beneath  the 
study  region.  These  are  preliminary  results,  so  we  invest  little  space  interpreting  the  variations.  We  started  the 
inversion  using  only  Rayleigh-wave  dispersion  in  the  period  range  from  7  to  200  seconds  (with  more  weight  on  the 
shorter  periods).  Figure  4  shows  example  fits  to  the  dispersion  values  and  the  Bouguer  grav  ity  variations.  As  seen  in 
earlier  studies  (Maeeira  and  Ammon,  2009)  the  addition  of  the  gravity  information  produces  little  change  to  the 
overall  fit  of  the  dispersion  values,  but  substantially  improves  the  fit  to  the  observed  surface  gravity  signals. 
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Observations 
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Figure  4.  Sample  dispersion  (top)  and  Bouguer  gravity  (bottom)  for  the  preliminary  inversion.  As  for  other 
regions  (Maeeira  and  Ammon,  2009),  the  inclusion  of  gravity  observations  produces  a  much  better 
Fit  to  the  observed  gravity  variations  with  little  change  to  the  Rayleigh-wave  Fit.  Note  that  the 
gravity  in  this  inversion  was  not  Filtered  to  emphasize  shorter  wavelength  variations  in  the  gravity 
to  focus  on  near-surface  variations. 


Software  Development  Progress 

A  large  component  of  the  proposed  w  ork  involves  improving  current  software  to  allow  more  easy  extension  to  larger 
focus  areas  and  more  flexible  incorporation  of  additional  data.  For  example,  the  earlier  work  in  Maeeira  and  Ammon 
(2009)  blended  gravity  and  surface-wave  observations  to  produce  a  shear-velocity/dcnsity  model  for  central  Asia, 
but  did  not  include  receiver  function  information.  Our  inversion  is  a  straightforw  ard  composite  of  coupled  and 
uncoupled  linearized  inversions  of  nonlinear  data  and  model  relationships.  In  essence  we  construct  a  large  set  of 
coupled  data  and  constrain  equations  into  the  form  Gm  =  d ,  where  the  model  parameters  are  in  m,  observations  and 
constraint  information  is  in  d ,  and  the  system  matrix,  G,  contains  partial  derivatives  and  coefficients  of  constraint 
equations. 

The  existing  inversion  tool  used  a  singular  value  decomposition  (SVD)  to  invert  the  system  matrix,  G.  The  SVD  is 
superb  numerical  choice  for  matrix  inversion  and  problem  analysis  but  the  simplest  implementations  require 
substantial  quantities  of  core  memory  and  computational  time  as  the  problem  grows  large.  To  reduce  these 
requirements  we  have  replaced  the  SVD  with  Page  and  Saunders  (1982)  LSQR  conjugate-gradient  based  routine. 

We  adopted  the  LSQR  function  from  the  Association  for  Computing  Machinery  (ACM)  FORTRAN  implementation 
last  modified  in  1994  (http://toms.acm.orK/).  Although  well  known  to  seismic  tomographers  (Pavlis,  1988;  Nolet, 
1993),  the  benefits  of  the  LSQR  algorithm  may  be  less  familiar  to  those  working  on  multi-dataset  inversion  tools. 
The  key  advantage  is  that  the  system  matrix,  G,  need  never  be  stored  in  memory.  G  can  be  stored  on  the  disk  as  it  is 
constructed  one  row  at  a  time.  The  LSQR  user  need  only  develop  a  subroutine  that  computes y  =  Gy  +  y  and 
x  =  Gtx  +  a*,  for  vectors  x  andy  provided  during  the  inversion  of  G.  G  can  be  written  to  disk  by  saving  only  the 
nonzero  columns  of  each  row  of  the  matrix  and  the  associated  indices  of  the  nonzero  columns.  Not  only  do  you  save 
storage  and  input-output  time,  you  also  avoid  waste  by  not  multiplying  numbers  by  zero.  In  our  early  experiments, 
we  have  saved  roughly  an  order  of  magnitude  in  computation  time  (measured  on  a  relatively  simple  desktop 
workstation).  The  great  reduction  in  computation  time  allows  more  experimentation  with  and  assessment  of  the 
model. 
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Example  Rayleigh-Wave  Dispersion  Receiver  Function  Coupled  Inversion 


Dispersion  Misfit  (km/s) 
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To  illustrate  the  result  of  combining 
surface-wave  dispersion  and 
receiver  functions  in  a  coupled 
inversion  of  cells  for  an 
approximation  to  a  3D  model,  we 
selected  a  small,  but  interesting 
region  of  the  focus  area  that 
includes  the  Caspian  Sea.  The 
Caspian  Sea  region  is  an  area  of 
complex  geologic  structure  (c.g., 
Mangino  and  Priestly,  1998)  and 
significantly  low  shear  velocities. 
The  sub  region  is  shown  in 
Figure  5.  which  also  shows  the 


Figure  5.  Initial  (left)  and  Final  (right)  dispersion  misfit  associated  with 
the  example  subregion  inversion  of  dispersion  only.  The  initial 
model  consists  of  a  simple  “continental'’  model  and  a  Caspian 
Sea  model.  The  initial  misfits  are  large,  but  the  convergence  is 
reasonable  to  mean  misfit  values  near  or  below  0.05  km/s. 


misfit  to  the  dispersion  values  for 
each  cell.  As  above,  the  dispersion 
values  are  from  the  University  of 
Colorado  tomographic  effort 
(Ritzwollcr  and  Levshin,  1998; 
Levshin  et  al.,  2001, 2002),  The 
initial  misfits  are  large,  but  the 

inversion  converges  nicely  to  reasonable  fits  after  four  iterations  (final)  Two  depth  slices  of  the  corresponding 
shear-velocity  model  (initial  and  final)  are  shown  in  Figure  6.  The  interpolation  and  contouring  algorithm  had  some 
trouble  with  the  rough  initial  model  that  contained  two  structures,  one  for  region  within  the  Caspian  Sea  and  another 
for  the  rest  of  the  model.  We  applied  a  simple  Laplacian  smoothing  betw  een  adjacent  cells  so  the  resulting  model 
remains  smooth. 


Initial,  Depth  =  5  km 


45'  50  55’ 


Initial,  Depth  -  37  km 


Figure  6.  Initial  (left)  and  final  (right)  shear-velocity  values 
associated  with  the  example  subregion  inversion  of  dispersion  only 
for  the  upper  and  lower  crust  or  upper  mantle  region  of  the  model. 
The  interpolation  distorts  the  rough  initial  model  producing 
artifacts  that  are  absent  when  the  inversion  is  complete  and  the 
model  is  smoother  (we  use  a  1°  cell  dimension  for  the  inversion). 


Next,  we  added  one  receiver  function  to 
the  dispersion  data  and  repeated  the  3D 
shear-velocity  inversion  The  receiver 
function  was  digitized  from  Mangino  and 
Priestly  ( 1 998),  and  includes  samples  the 
structure  in  our  cell  with  the  lower  left 
comer  at  53E  and  40N  (Figure  7). 

Figure  7  is  a  plot  of  the  initial  and  final 
dispersion  fits,  which  show  little  change 
from  those  without  the  receiver  function. 
The  initial  and  final  fit  to  the  receiver 
function  are  shown  in  Figure  8.  The  effect 
of  including  the  receiver  function  is  to 
thicken  the  erust  and  better  define  the 
crust-mantle  boundary.  The  net  result  is  a 
decrease  in  the  velocity  at  the  depth  slice 
near  37  km.  Two  depth  slices  of  the 
corresponding  shear-velocity  model 
(initial  and  final)  arc  shown  in  Figure  9. 
The  thickening  of  the  crust  near  the 
receiver  function  observation  is  clearly 
visible  in  the  deeper  slice  through  the 
model  (compare  Figure  6  and  Figure  9). 
The  fit  to  the  receiver  function  is  not 
perfect  and  the  spread  of  the  slower  deep 
crustal  speeds  into  the  Caspian  region 
docs  not  make  sense  geologically. 
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These  effects  are  driven  by  the  simple  L  aplaeian  smoothness  requirements  placed  on  the  3D  model. 
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Figure  7.  Initial  (left)  and  final  (right)  dispersion  misfit  associated  with  the  example  subregion  inversion  of 
dispersion  and  a  receiver  function.  The  initial  model  consists  of  a  simple  “continental”  model  and  a 
Caspian  Sea  model.  The  initial  misfits  are  large,  but  the  convergence  is  reasonable  to  mean  misfit 
values  near  or  below  0.05  km/s.  The  cell  containing  the  isolated  receiver  function  is  labeled. 


Figure  8.  Observed  (top),  final  prediction  (middle)  and  initial  prediction  (bottom)  of  the  receiver  functions 

from  the  inversion  shown  in  Figure  4.  Vertical  lines  show  the  arrivals  most  likclv  associated  with  the 
crust-mantle  transition.  Although  the  final  fit  is  certainly  not  perfect,  note  the  improvement  in  timing 
of  the  main  crustal  conv  ersion  and  rev  erberation  (v  ertical  lines)  from  the  initial  to  final  models. 

Including  Geologic  Information 

To  produce  models  that  have  realistic  ‘sharp’  boundaries  requires  that  we  include  independent  information  on  the 
location  of  those  boundaries.  Such  information  is  available  (for  the  shallow  part  of  the  model)  in  independent  data 
sets  such  as  gravity,  surface  geologic  maps,  and  even  something  as  simple  as  topography.  As  part  of  this  w  ork,  we 
plan  to  resolve  sharp  features  by  adapting  our  imaging  algorithms  to  allow  the  inclusion  of  geologic  information  on 
the  location  and  nature  of  the  boundary  into  shear-velocity  inversions  that  permit  such  features  (implemented 
through  custom  geologic  smoothness  constraints  that  allow  velocities  to  be  de-eorrelated  across  major  geologic 
transitions).  Including  a  priori  information  into  an  inversion  is  obviously  only  as  good  as  the  information  that  is 
included.  Thus  the  inclusion  of  this  type  of  information  into  the  reconstruction  of  shear- velocity  models  of  the 
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subsurface  must  proceed  carefully  and  include  documentation  of  the  importance  of  the  assumed  a  priori  information 
on  the  resulting  model. 


Final.  Depth  =  5  km 


Shear  Velocity  (knVij 


3.5 

3.0 

2.5 


Figure  9.  Initial  (left)  and  final  (right)  interpolated  shear-*  eloeit\  values  associated  with  the  example 

subregion  inversion  of  dispersion  and  a  receiver  function.  The  addition  of  the  receiver  function 
improves  the  estimate  of  crustal  details.  Please  note  the  dramatic  difference  in  the  \eloeity  scales 
shown  to  the  right  of  each  pair  of  images.  The  interpolation  distorts  the  rough  initial  model 
producing  artifacts  that  arc  absent  w  hen  the  imersion  is  complete  and  the  model  is  smoother 
(we  use  a  1°  cell  dimension  for  the  inversion). 

We  have  made  progress  on  including  this  type  of  information  into  the  inversion,  but  have  not  finished  debugging 
this  component  of  the  softw  are.  We  expect  that  uncoupling  the  structures  beneath  a  region  of  the  Caspian  Sea  and  its 
eastern  margin  will  for  example,  allow  a  better  fit  to  the  KRF  receiver  function  and  decrease  the  leaking  of  receiver 
function  model  information  into  the  model  beneath  the  Caspian  Sea.  Implementation  of  these  smoothing  constraints 
requires  more  work,  but  in  the  end  includes  more  information  in  the  shear- velocity  model  reconstruction,  and 
hopefully  leads  to  improved  regional  earth  models. 

CONCLUSIONS 


We  have  initiated  a  two-year  project  to  map  the  subsurface  geologic  variations  using  seismic  dispersion,  gravity,  and 
receiver-function  observations.  We  face  significant  challenges  in  our  efforts  to  include  effective  point  constraints  on 
structure  (receiver  functions)  with  the  spatially  continuous  surface-wave  tomography  and  gravity  observations.  Our 
work  complements  ongoing  work  at  LANL  to  integrate  body-wave  travel  times  into  the  same  formalism.  The  basic 
philosophy  is  that  models  that  explain  more  data  are  better.  The  ultimate  utility  of  the  derived  earth  models  is  to 
provide  improved  predictive  capabilities  for  routine  seismic  analyses  and  to  provide  adequate  starting  models  for  3D 
waveform  inversion  approaches. 
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